1.Distances et dissimilarites

Sokal et Michener

SokMichener = function(i,j){ #i et j deux individus (vecteurs)  
  p    = length(i)
  diss = sum(i!=j)/p
  return(diss)
}

Donnees:

ExBinaire = read.table("donnees/exemple.dat") #27 individus et 19 variables binaires

MatDiss = outer(as.data.frame(t(ExBinaire)),as.data.frame(t(ExBinaire)),Vectorize(SokMichener))
MatDiss = as.dist(MatDiss)

Distance chi-deux

#install.packages("ade4")
library(ade4)
## Warning: package 'ade4' was built under R version 3.4.4
accidents = read.table("donnees/accidents.txt",header=TRUE)
acc       = accidents[1001:2000,]
p         = ncol(acc)

Transformation des variables Nb_Blegers et Nb_BGrave en variables binaires

acc$Nb_Blegers = as.numeric(acc$Nb_Blegers!=0) 
acc$Nb_BGrave = as.numeric(acc$Nb_BGrave!=0)

Recherche de variables manquantes

(0 dans les 4 premieres colonnes)

for (i in 1:4) acc[acc[,i]==0,i]=NA #concerne uniquement la variable n°3
acc = na.omit(acc) #on retire les individus avec des donnees manquantes

acc.disj = as.matrix(acm.disjonctif(acc)) #transformation en matrice pour accelerer les calculs

Fonction

chi2.dist = function(acc.disj,i,j){ #acc.disj : tableau disjonctif; i, j : indice des individus
  njl = apply(acc.disj,2,sum) # effectif modalites
  delta = abs(acc.disj[i,]-acc.disj[j,])
  dist = sum(delta/njl) #distance du chi-deux a une constante multiplicative pres
  return(dist)
}

Optimiser le code

Pour optimiser le temps de calcul, on n’utilise pas la fonction definie au (4).

t = proc.time()
njl = apply(acc.disj,2,sum)
distk2 = sapply(1:nrow(acc),FUN=function(j){
  apply(sapply(1:nrow(acc),FUN=function(i) abs(acc.disj[j,]-acc.disj[i,]),simplify = TRUE)/njl,2,sum)
})
proc.time()-t ## 14 secondes
##    user  system elapsed 
##    5.93    0.00    5.96
distk2 = as.dist(distk2)

Distance de correlation

Distance utile lorsqu’on s’interesse aux proximites de formes plutot que de valeurs:

LoupChien = read.table("donnees/loup_chien.dat",header=TRUE)
LCdata    = LoupChien[LoupChien$espece!="inconnu",1:3]
LCespece  = LoupChien$espece[LoupChien$espece!="inconnu"]

head(LCdata)
##         LoC  LoM  LaM
## loup1  23.4 15.0 19.1
## chien1 17.3 12.8 14.3
## chien2 20.5 13.7 16.6
## loup2  26.0 16.0 19.4
## chien3 20.5 14.4 17.7
## chien4 19.0 12.4 14.9
CorDist = as.dist(cor(t(LCdata)))
head(as.matrix(CorDist))
##            loup1    chien1    chien2     loup2    chien3    chien4
## loup1  0.0000000 0.9844853 0.9974837 0.9857625 0.9981380 0.9921584
## chien1 0.9844853 0.0000000 0.9944480 0.9999724 0.9719494 0.9986965
## chien2 0.9974837 0.9944480 0.0000000 0.9952027 0.9913020 0.9985229
## loup2  0.9857625 0.9999724 0.9952027 0.0000000 0.9736708 0.9990483
## chien3 0.9981380 0.9719494 0.9913020 0.9736708 0.0000000 0.9826872
## chien4 0.9921584 0.9986965 0.9985229 0.9990483 0.9826872 0.0000000
##           chien5    chien6    chien7     loup3     loup4    chien8
## loup1  0.9223176 0.9110901 0.9472144 0.9988106 0.9993217 0.9400963
## chien1 0.9758143 0.9692835 0.9887735 0.9918698 0.9773556 0.9853293
## chien2 0.9473934 0.9380214 0.9675603 0.9997541 0.9941961 0.9618999
## loup2  0.9741624 0.9674285 0.9876355 0.9927883 0.9789015 0.9840334
## chien3 0.8970293 0.8842505 0.9258953 0.9939768 0.9997072 0.9175517
## chien4 0.9633842 0.9554662 0.9798576 0.9970724 0.9868824 0.9753336
##           chien9   chien10     loup5   chien11   chien12   chien13
## loup1  0.9462001 0.9854723 0.7495379 0.9752034 0.9621903 0.9832804
## chien1 0.9882981 0.9999837 0.6217568 0.9989061 0.9950556 0.9999773
## chien2 0.9667599 0.9950333 0.7007213 0.9884396 0.9790797 0.9937162
## loup2  0.9871369 0.9999985 0.6275616 0.9985309 0.9942898 0.9998996
## chien3 0.9247010 0.9732779 0.7885194 0.9598885 0.9437847 0.9703423
## chien4 0.9792239 0.9989719 0.6609237 0.9952172 0.9886889 0.9983298
##            loup6   chien14   chien15   chien16   chien17   chien18
## loup1  0.9986409 0.9980399 0.9944776 0.9863727 0.9999055 0.9962420
## chien1 0.9922924 0.9935365 0.9974637 0.9999384 0.9819805 0.9959834
## chien2 0.9998230 0.9999653 0.9994157 0.9955550 0.9964150 0.9998757
## loup2  0.9931861 0.9943529 0.9979652 0.9999933 0.9833582 0.9966215
## chien3 0.9936024 0.9923643 0.9862244 0.9745006 0.9988821 0.9891039
## chien4 0.9973241 0.9980354 0.9997966 0.9992016 0.9903467 0.9992554
##          chien19   chien20     loup7   chien21     loup8     loup9
## loup1  0.9426886 0.9902794 0.9767248 0.9966248 0.8066823 0.9894194
## chien1 0.9866118 0.9993217 0.9992082 0.9955668 0.8978653 0.9995263
## chien2 0.9639727 0.9976487 0.9894740 0.9999370 0.8465511 0.9972156
## loup2  0.9853722 0.9995678 0.9988849 0.9962385 0.8945677 0.9997274
## chien3 0.9205804 0.9799515 0.9618226 0.9897618 0.7691324 0.9787275
## chien4 0.9770013 0.9998988 0.9958749 0.9990700 0.8742221 0.9997943
##          chien22   chien23    loup10    loup11   chien24   chien25
## loup1  0.9999908 0.9953919 0.9838278 0.9858011 0.9702432 0.9976516
## chien1 0.9852298 0.9631231 0.9999931 0.9999706 0.9976764 0.9941916
## chien2 0.9977790 0.9860889 0.9940509 0.9952252 0.9849681 0.9999971
## loup2  0.9864755 0.9650966 0.9999379 1.0000000 0.9971424 0.9949642
## chien3 0.9978668 0.9993874 0.9710704 0.9737233 0.9536675 0.9916161
## chien4 0.9926860 0.9756014 0.9985003 0.9990583 0.9928984 0.9983891
##          chien26   chien27    loup12   chien28   chien29   chien30
## loup1  0.9605669 0.9843357 0.9666393 0.9997508 0.9967406 0.9773556
## chien1 0.9944524 0.9381285 0.9965866 0.9881567 0.9954320 0.9993217
## chien2 0.9778625 0.9693595 0.9823664 0.9988177 0.9999519 0.9898981
## loup2  0.9936430 0.9406767 0.9959454 0.9892700 0.9961142 0.9990203
## chien3 0.9418184 0.9932568 0.9492157 0.9965278 0.9899639 0.9626287
## chien4 0.9877870 0.9545812 0.9910736 0.9947011 0.9990076 0.9961392

2.Classifiation Automatique

Donnees Iris

summary(iris)
##   Sepal.Length    Sepal.Width     Petal.Length    Petal.Width   
##  Min.   :4.300   Min.   :2.000   Min.   :1.000   Min.   :0.100  
##  1st Qu.:5.100   1st Qu.:2.800   1st Qu.:1.600   1st Qu.:0.300  
##  Median :5.800   Median :3.000   Median :4.350   Median :1.300  
##  Mean   :5.843   Mean   :3.057   Mean   :3.758   Mean   :1.199  
##  3rd Qu.:6.400   3rd Qu.:3.300   3rd Qu.:5.100   3rd Qu.:1.800  
##  Max.   :7.900   Max.   :4.400   Max.   :6.900   Max.   :2.500  
##        Species  
##  setosa    :50  
##  versicolor:50  
##  virginica :50  
##                 
##                 
## 
help(iris)
## starting httpd help server ... done
plot(iris)

summary(iris[1:50,])
##   Sepal.Length    Sepal.Width     Petal.Length    Petal.Width   
##  Min.   :4.300   Min.   :2.300   Min.   :1.000   Min.   :0.100  
##  1st Qu.:4.800   1st Qu.:3.200   1st Qu.:1.400   1st Qu.:0.200  
##  Median :5.000   Median :3.400   Median :1.500   Median :0.200  
##  Mean   :5.006   Mean   :3.428   Mean   :1.462   Mean   :0.246  
##  3rd Qu.:5.200   3rd Qu.:3.675   3rd Qu.:1.575   3rd Qu.:0.300  
##  Max.   :5.800   Max.   :4.400   Max.   :1.900   Max.   :0.600  
##        Species  
##  setosa    :50  
##  versicolor: 0  
##  virginica : 0  
##                 
##                 
## 
summary(iris[51:100,])
##   Sepal.Length    Sepal.Width     Petal.Length   Petal.Width   
##  Min.   :4.900   Min.   :2.000   Min.   :3.00   Min.   :1.000  
##  1st Qu.:5.600   1st Qu.:2.525   1st Qu.:4.00   1st Qu.:1.200  
##  Median :5.900   Median :2.800   Median :4.35   Median :1.300  
##  Mean   :5.936   Mean   :2.770   Mean   :4.26   Mean   :1.326  
##  3rd Qu.:6.300   3rd Qu.:3.000   3rd Qu.:4.60   3rd Qu.:1.500  
##  Max.   :7.000   Max.   :3.400   Max.   :5.10   Max.   :1.800  
##        Species  
##  setosa    : 0  
##  versicolor:50  
##  virginica : 0  
##                 
##                 
## 
summary(iris[101:150,])
##   Sepal.Length    Sepal.Width     Petal.Length    Petal.Width   
##  Min.   :4.900   Min.   :2.200   Min.   :4.500   Min.   :1.400  
##  1st Qu.:6.225   1st Qu.:2.800   1st Qu.:5.100   1st Qu.:1.800  
##  Median :6.500   Median :3.000   Median :5.550   Median :2.000  
##  Mean   :6.588   Mean   :2.974   Mean   :5.552   Mean   :2.026  
##  3rd Qu.:6.900   3rd Qu.:3.175   3rd Qu.:5.875   3rd Qu.:2.300  
##  Max.   :7.900   Max.   :3.800   Max.   :6.900   Max.   :2.500  
##        Species  
##  setosa    : 0  
##  versicolor: 0  
##  virginica :50  
##                 
##                 
##

CAH

dist.Iris <- dist(iris[,1:4])
cah.ward.Iris = hclust(dist.Iris,method="ward.D2")
plot(cah.ward.Iris,hang=-1)

Choix du nmb de groupes

plot(rev(cah.ward.Iris$height),type="b", main="Courbe : perte de l'inertie")
# Choix du nombre de groupes de 3 ou 4 (moment du premier "grand saut")
abline(h=c(3,4), col=c("blue","red"))

gpe.ward.Iris3 = cutree(cah.ward.Iris,k=3)
gpe.ward.Iris4 = cutree(cah.ward.Iris,k=4)
plot(iris, pch=gpe.ward.Iris3, col=gpe.ward.Iris3)

plot(iris, pch=gpe.ward.Iris4, col=gpe.ward.Iris4)

groupes3=gpe.ward.Iris3
groupes4=gpe.ward.Iris4
#on va choisir 3 groupes car cela à l'air mieux rpz même s'il y a que une espèce très bien rpz et les deux autres se moins biens séparé

Groupes: 3 & 4

plot(cah.ward.Iris,hang=-1)
rect.hclust(cah.ward.Iris, 3, border ="blue")

plot(cah.ward.Iris,hang=-1)
rect.hclust(cah.ward.Iris, 4, border ="blue")

table(iris$Species, groupes3)
##             groupes3
##               1  2  3
##   setosa     50  0  0
##   versicolor  0 49  1
##   virginica   0 15 35
#Après CAH, l'espèce setosa forme un groupe entier à lui même contrairement au deux autres espèces (même si chaque groupe est majoritairement composé d'une espèce)
table(iris$Species, groupes4)
##             groupes4
##               1  2  3  4
##   setosa     50  0  0  0
##   versicolor  0 24 25  1
##   virginica   0 14  1 35
#De même pour la l'espèce setosa, les autres groupes sont formés des deux autres espèce mais cette fois sans avec un groupe de majorité 'extrême' présente par une espèces (groupe 2)

Cusplots

library(cluster)
## Warning: package 'cluster' was built under R version 3.4.4
help("clusplot")
clusplot(iris[,1:4],groupes3, main="Clusplot Iris")

clusplot(iris[,1:4],groupes4, main="Clusplot Iris")

#cette rpz correspond au graphique de partition des Iris (la répartition par groupes), et la distance entre ses groupes.

3.Agrégation autour de centres mobiles

K-means

#help(kmeans)
#paramètres : la bd, le nmb de groupes, le nmb d'itératin max de l'algorithme, nstars grand si données très grandes : ici 150 donc 10 au max c'est bien
IrisKmeans <- kmeans(x=iris[,1:4],centers = 3, nstart = 50)
Kgroupes <- as.factor(IrisKmeans$cluster)
IrisKmeans <- cbind(iris, Kgroupes)
IrisKmeans <- as.data.frame(IrisKmeans)
library(plotly)
plot_ly(data=IrisKmeans,  
        x=IrisKmeans$Sepal.Length, 
        y=IrisKmeans$Sepal.Width, 
        z=IrisKmeans$Petal.Length, 
        color=IrisKmeans$Kgroupes, 
        colors = c("royalblue4", "gold", "darkgrey"))
library(ggplot2)
ggplot(IrisKmeans)+aes(x=IrisKmeans$Sepal.Length, fill=IrisKmeans$Kgroupes)+geom_bar()+labs(x="Sepal Length", y="Effectif", fill="Groupes")

ggplot(IrisKmeans)+aes(x=IrisKmeans$Petal.Length, fill=IrisKmeans$Kgroupes)+geom_bar()+labs(x="Pepal Length", y="Effectif", fill="Groupes")

ggplot(IrisKmeans)+aes(x=IrisKmeans$Sepal.Width, fill=IrisKmeans$Kgroupes)+geom_bar()+labs(x="Sepal Width", y="Effectif", fill="Groupes")

ggplot(IrisKmeans)+aes(x=IrisKmeans$Petal.Width, fill=IrisKmeans$Kgroupes)+geom_bar()+labs(x="Pepal Width", y="Effectif", fill="Groupes")

#le dernier graphique rpz bien la séparation du groupe 1 avec les autres groupes
#le graph 2 et 4 rpz la séparation entre les deux autres groupes

#en relancant plusieurs fois la fonction on remarque que les noms des groupes changent : cad que l'algorithme définit de facon aléatoire les groupes mais leurs caractéristiques sont les mêmes (plus au moins) : au final on gardera la meilleur classification : et c'est grace à l'option 'nstart' (nomnre d'initialisation aléatoire à essayer) c'est mieux d'en utiliser plus